

<!DOCTYPE html>
<!--[if IE 8]><html class="no-js lt-ie9" lang="en" > <![endif]-->
<!--[if gt IE 8]><!--> <html class="no-js" lang="en" > <!--<![endif]-->
<head>
  <meta charset="utf-8">
  
  <meta name="viewport" content="width=device-width, initial-scale=1.0">
  
  <title>quippy.structures &mdash; quippy 7ea04c5+ documentation</title>
  

  
  
    <link rel="shortcut icon" href="../../_static/favicon.ico"/>
  
  
  

  
  <script type="text/javascript" src="../../_static/js/modernizr.min.js"></script>
  
    
      <script type="text/javascript" id="documentation_options" data-url_root="../../" src="../../_static/documentation_options.js"></script>
        <script type="text/javascript" src="../../_static/jquery.js"></script>
        <script type="text/javascript" src="../../_static/underscore.js"></script>
        <script type="text/javascript" src="../../_static/doctools.js"></script>
        <script type="text/javascript" src="../../_static/language_data.js"></script>
        <script async="async" type="text/javascript" src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.5/latest.js?config=TeX-AMS-MML_HTMLorMML"></script>
        <script type="text/x-mathjax-config">MathJax.Hub.Config({"tex2jax": {"processClass": "math|output_area", "processEscapes": true, "ignoreClass": "document", "inlineMath": [["$", "$"], ["\\(", "\\)"]]}})</script>
    
    <script type="text/javascript" src="../../_static/js/theme.js"></script>

    

  
  <link rel="stylesheet" href="../../_static/css/theme.css" type="text/css" />
  <link rel="stylesheet" href="../../_static/pygments.css" type="text/css" />
  <link rel="stylesheet" href="../../_static/theme_overrides.css" type="text/css" />
    <link rel="index" title="Index" href="../../genindex.html" />
    <link rel="search" title="Search" href="../../search.html" /> 
</head>

<body class="wy-body-for-nav">

   
  <div class="wy-grid-for-nav">
    
    <nav data-toggle="wy-nav-shift" class="wy-nav-side">
      <div class="wy-side-scroll">
        <div class="wy-side-nav-search" >
          

          
            <a href="../../index.html" class="icon icon-home"> quippy
          

          
            
            <img src="../../_static/hybrid.png" class="logo" alt="Logo"/>
          
          </a>

          
            
            
              <div class="version">
                7ea04c5+
              </div>
            
          

          
<div role="search">
  <form id="rtd-search-form" class="wy-form" action="../../search.html" method="get">
    <input type="text" name="q" placeholder="Search docs" />
    <input type="hidden" name="check_keywords" value="yes" />
    <input type="hidden" name="area" value="default" />
  </form>
</div>

          
        </div>

        <div class="wy-menu wy-menu-vertical" data-spy="affix" role="navigation" aria-label="main navigation">
          
    
            
            
              
            
            
              <ul>
<li class="toctree-l1"><a class="reference internal" href="../../intro.html">Introduction to <code class="docutils literal notranslate"><span class="pre">quippy</span></code></a></li>
<li class="toctree-l1"><a class="reference internal" href="../../install.html">Installation of QUIP and quippy</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../Tutorials/index.html">Tutorials</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../quippy.html">Quippy library reference</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../visualisation.html">Visualisation Tools</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../fortran_wrapper.html">Appendix: wrapping Fortran 90 Code</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../teach_sparse.html">Command line options of the teach_sparse main program</a></li>
</ul>

            
          
    <p></p>
    <ul>
        <li><a href="../../genindex.html">Index</a></li>
        <li><a href="../../py-modindex.html">Module Index</a></li>
        <li><a href="../../search.html">Search Page</a></li>
        <li><a href="../index.html">Module Source Code</a></li>
    </ul>
  
        </div>
      </div>
    </nav>

    <section data-toggle="wy-nav-shift" class="wy-nav-content-wrap">

      
      <nav class="wy-nav-top" aria-label="top navigation">
        
          <i data-toggle="wy-nav-top" class="fa fa-bars"></i>
          <a href="../../index.html">quippy</a>
        
      </nav>


      <div class="wy-nav-content">
        
        <div class="rst-content">
        
          















<div role="navigation" aria-label="breadcrumbs navigation">

  <ul class="wy-breadcrumbs">
    
      <li><a href="../../index.html">Docs</a> &raquo;</li>
        
          <li><a href="../index.html">Module code</a> &raquo;</li>
        
      <li>quippy.structures</li>
    
    
      <li class="wy-breadcrumbs-aside">
        
      </li>
    
  </ul>

  
  <hr/>
</div>
          <div role="main" class="document" itemscope="itemscope" itemtype="http://schema.org/Article">
           <div itemprop="articleBody">
            
  <h1>Source code for quippy.structures</h1><div class="highlight"><pre>
<span></span><span class="c1"># HQ XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX</span>
<span class="c1"># HQ X</span>
<span class="c1"># HQ X   quippy: Python interface to QUIP atomistic simulation library</span>
<span class="c1"># HQ X</span>
<span class="c1"># HQ X   Copyright James Kermode 2010</span>
<span class="c1"># HQ X</span>
<span class="c1"># HQ X   These portions of the source code are released under the GNU General</span>
<span class="c1"># HQ X   Public License, version 2, http://www.gnu.org/copyleft/gpl.html</span>
<span class="c1"># HQ X</span>
<span class="c1"># HQ X   If you would like to license the source code under different terms,</span>
<span class="c1"># HQ X   please contact James Kermode, james.kermode@gmail.com</span>
<span class="c1"># HQ X</span>
<span class="c1"># HQ X   When using this software, please cite the following reference:</span>
<span class="c1"># HQ X</span>
<span class="c1"># HQ X   http://www.jrkermode.co.uk/quippy</span>
<span class="c1"># HQ X</span>
<span class="c1"># HQ XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX</span>

<span class="sd">&quot;&quot;&quot;</span>
<span class="sd">This module contains a variety of structure generating routines.</span>
<span class="sd">&quot;&quot;&quot;</span>

<span class="kn">import</span> <span class="nn">numpy</span> <span class="k">as</span> <span class="nn">np</span>
<span class="kn">import</span> <span class="nn">itertools</span>
<span class="kn">from</span> <span class="nn">quippy.atoms</span> <span class="k">import</span> <span class="o">*</span>
<span class="kn">import</span> <span class="nn">quippy._structures</span>
<span class="kn">from</span> <span class="nn">quippy._structures</span> <span class="k">import</span> <span class="o">*</span>
<span class="kn">from</span> <span class="nn">quippy</span> <span class="k">import</span> <span class="n">quippy_array</span>
<span class="kn">from</span> <span class="nn">quippy.farray</span> <span class="k">import</span> <span class="n">FortranArray</span><span class="p">,</span> <span class="n">fidentity</span><span class="p">,</span> <span class="n">fzeros</span><span class="p">,</span> <span class="n">frange</span><span class="p">,</span> <span class="n">farray</span><span class="p">,</span> <span class="n">gcd</span>

<span class="kn">import</span> <span class="nn">numpy</span> <span class="k">as</span> <span class="nn">np</span>

<span class="n">__all__</span> <span class="o">=</span> <span class="n">quippy</span><span class="o">.</span><span class="n">_structures</span><span class="o">.</span><span class="n">__all__</span> <span class="o">+</span> <span class="p">[</span><span class="s1">&#39;orthorhombic_slab&#39;</span><span class="p">,</span> <span class="s1">&#39;rotation_matrix&#39;</span><span class="p">,</span>
                                        <span class="s1">&#39;quartz_params&#39;</span><span class="p">,</span> <span class="s1">&#39;get_bulk_params&#39;</span><span class="p">,</span>
                                        <span class="s1">&#39;get_quartz_params&#39;</span><span class="p">,</span> <span class="s1">&#39;get_bond_lengths&#39;</span><span class="p">,</span>
                                        <span class="s1">&#39;MillerIndex&#39;</span><span class="p">,</span> <span class="s1">&#39;MillerPlane&#39;</span><span class="p">,</span> <span class="s1">&#39;MillerDirection&#39;</span><span class="p">,</span> <span class="s1">&#39;angle_between&#39;</span><span class="p">]</span>

<div class="viewcode-block" id="MillerIndex"><a class="viewcode-back" href="../../structures.html#quippy.structures.MillerIndex">[docs]</a><span class="k">class</span> <span class="nc">MillerIndex</span><span class="p">(</span><span class="n">quippy_array</span><span class="p">):</span>
    <span class="sd">&quot;&quot;&quot;</span>
<span class="sd">    Representation of a three of four index Miller direction or plane</span>

<span class="sd">    A :class:`MillerIndex` can be constructed from vector or parsed from a string::</span>

<span class="sd">        x = MillerIndex(&#39;-211&#39;)</span>
<span class="sd">        y = MillerIndex(&#39;111&#39;, type=&#39;plane&#39;)</span>
<span class="sd">        z = x.cross(y)</span>
<span class="sd">        print x # prints &quot;[-211]&quot;</span>
<span class="sd">        print y # prints &quot;(111)&quot;, note round brackets denoting a plane</span>
<span class="sd">        print z.latex()</span>
<span class="sd">        assert(angle_between(x,y) == pi/2.)</span>
<span class="sd">        assert(angle_between(y,z) == pi/2.)</span>
<span class="sd">        assert(angle_between(x,z) == pi/2.)</span>
<span class="sd">    &quot;&quot;&quot;</span>

    <span class="n">__array_priority__</span> <span class="o">=</span> <span class="mf">101.0</span>

    <span class="n">brackets</span> <span class="o">=</span> <span class="p">{</span><span class="s1">&#39;direction&#39;</span><span class="p">:</span> <span class="s1">&#39;[]&#39;</span><span class="p">,</span>
                <span class="s1">&#39;direction_family&#39;</span><span class="p">:</span> <span class="s1">&#39;&lt;&gt;&#39;</span><span class="p">,</span>
                <span class="s1">&#39;plane&#39;</span><span class="p">:</span> <span class="s1">&#39;()&#39;</span><span class="p">,</span>
                <span class="s1">&#39;plane_family&#39;</span><span class="p">:</span> <span class="s1">&#39;</span><span class="si">{}</span><span class="s1">&#39;</span><span class="p">}</span>

    <span class="n">all_brackets</span> <span class="o">=</span> <span class="nb">list</span><span class="p">(</span><span class="n">itertools</span><span class="o">.</span><span class="n">chain</span><span class="p">(</span><span class="o">*</span><span class="n">brackets</span><span class="o">.</span><span class="n">values</span><span class="p">()))</span>

    <span class="k">def</span> <span class="nf">__new__</span><span class="p">(</span><span class="bp">cls</span><span class="p">,</span> <span class="n">v</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span> <span class="nb">type</span><span class="o">=</span><span class="s1">&#39;direction&#39;</span><span class="p">):</span>
        <span class="k">if</span> <span class="nb">isinstance</span><span class="p">(</span><span class="n">v</span><span class="p">,</span> <span class="n">basestring</span><span class="p">):</span>
            <span class="n">v</span> <span class="o">=</span> <span class="n">MillerIndex</span><span class="o">.</span><span class="n">parse</span><span class="p">(</span><span class="n">v</span><span class="p">)</span>
        <span class="k">if</span> <span class="nb">len</span><span class="p">(</span><span class="n">v</span><span class="p">)</span> <span class="o">==</span> <span class="mi">3</span> <span class="ow">or</span> <span class="nb">len</span><span class="p">(</span><span class="n">v</span><span class="p">)</span> <span class="o">==</span> <span class="mi">4</span><span class="p">:</span>
            <span class="bp">self</span> <span class="o">=</span> <span class="n">quippy_array</span><span class="o">.</span><span class="fm">__new__</span><span class="p">(</span><span class="bp">cls</span><span class="p">,</span> <span class="n">v</span><span class="p">)</span>
        <span class="k">else</span><span class="p">:</span>
            <span class="k">raise</span> <span class="ne">ValueError</span><span class="p">(</span><span class="s1">&#39;</span><span class="si">%s</span><span class="s1"> input v should be of length 3 or 4&#39;</span> <span class="o">%</span> <span class="bp">cls</span><span class="o">.</span><span class="vm">__name__</span><span class="p">)</span>
        <span class="bp">self</span><span class="o">.</span><span class="n">type</span> <span class="o">=</span> <span class="nb">type</span>
        <span class="bp">self</span><span class="o">.</span><span class="n">simplify</span><span class="p">()</span>
        <span class="k">return</span> <span class="bp">self</span>

    <span class="k">def</span> <span class="nf">__array_finalize__</span><span class="p">(</span><span class="bp">self</span><span class="p">,</span> <span class="n">obj</span><span class="p">):</span>
        <span class="k">if</span> <span class="n">obj</span> <span class="ow">is</span> <span class="kc">None</span><span class="p">:</span>
            <span class="k">return</span>
        <span class="bp">self</span><span class="o">.</span><span class="n">type</span> <span class="o">=</span> <span class="nb">getattr</span><span class="p">(</span><span class="n">obj</span><span class="p">,</span> <span class="s1">&#39;type&#39;</span><span class="p">,</span> <span class="s1">&#39;direction&#39;</span><span class="p">)</span>

    <span class="k">def</span> <span class="nf">__repr__</span><span class="p">(</span><span class="bp">self</span><span class="p">):</span>
        <span class="k">return</span> <span class="p">(</span><span class="s1">&#39;</span><span class="si">%s</span><span class="s1">([&#39;</span><span class="o">+</span><span class="s1">&#39;</span><span class="si">%d</span><span class="s1">&#39;</span><span class="o">*</span><span class="nb">len</span><span class="p">(</span><span class="bp">self</span><span class="p">)</span><span class="o">+</span><span class="s1">&#39;])&#39;</span><span class="p">)</span> <span class="o">%</span> <span class="p">((</span><span class="bp">self</span><span class="o">.</span><span class="vm">__class__</span><span class="o">.</span><span class="vm">__name__</span><span class="p">,)</span> <span class="o">+</span> <span class="nb">tuple</span><span class="p">(</span><span class="bp">self</span><span class="p">))</span>
        
    <span class="k">def</span> <span class="nf">__str__</span><span class="p">(</span><span class="bp">self</span><span class="p">):</span>
        <span class="n">bopen</span><span class="p">,</span> <span class="n">bclose</span> <span class="o">=</span> <span class="n">MillerIndex</span><span class="o">.</span><span class="n">brackets</span><span class="p">[</span><span class="bp">self</span><span class="o">.</span><span class="n">type</span><span class="p">]</span>
        <span class="k">return</span> <span class="p">(</span><span class="n">bopen</span><span class="o">+</span><span class="s1">&#39;</span><span class="si">%d</span><span class="s1">&#39;</span><span class="o">*</span><span class="nb">len</span><span class="p">(</span><span class="bp">self</span><span class="p">)</span><span class="o">+</span><span class="n">bclose</span><span class="p">)</span> <span class="o">%</span> <span class="nb">tuple</span><span class="p">(</span><span class="bp">self</span><span class="p">)</span>

<div class="viewcode-block" id="MillerIndex.latex"><a class="viewcode-back" href="../../structures.html#quippy.structures.MillerIndex.latex">[docs]</a>    <span class="k">def</span> <span class="nf">latex</span><span class="p">(</span><span class="bp">self</span><span class="p">):</span>
        <span class="sd">&quot;&quot;&quot;</span>
<span class="sd">        Format this :class:`MillerIndex` as a LaTeX string</span>
<span class="sd">        &quot;&quot;&quot;</span>
        <span class="n">s</span> <span class="o">=</span> <span class="s1">&#39;$&#39;</span>
        <span class="n">bopen</span><span class="p">,</span> <span class="n">bclose</span> <span class="o">=</span> <span class="n">MillerIndex</span><span class="o">.</span><span class="n">brackets</span><span class="p">[</span><span class="bp">self</span><span class="o">.</span><span class="n">type</span><span class="p">]</span>
        <span class="n">s</span> <span class="o">+=</span> <span class="n">bopen</span>
        <span class="k">for</span> <span class="n">component</span> <span class="ow">in</span> <span class="bp">self</span><span class="p">:</span>
            <span class="k">if</span> <span class="n">component</span> <span class="o">&lt;</span> <span class="mi">0</span><span class="p">:</span>
                <span class="n">s</span> <span class="o">+=</span> <span class="sa">r</span><span class="s1">&#39;\bar{</span><span class="si">%d</span><span class="s1">}&#39;</span> <span class="o">%</span> <span class="nb">abs</span><span class="p">(</span><span class="n">component</span><span class="p">)</span>
            <span class="k">else</span><span class="p">:</span>
                <span class="n">s</span> <span class="o">+=</span> <span class="s1">&#39;</span><span class="si">%d</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="n">component</span>
        <span class="n">s</span> <span class="o">+=</span> <span class="n">bclose</span>
        <span class="n">s</span> <span class="o">+=</span> <span class="s1">&#39;$&#39;</span>
        <span class="k">return</span> <span class="n">s</span></div>

<div class="viewcode-block" id="MillerIndex.parse"><a class="viewcode-back" href="../../structures.html#quippy.structures.MillerIndex.parse">[docs]</a>    <span class="nd">@classmethod</span>
    <span class="k">def</span> <span class="nf">parse</span><span class="p">(</span><span class="bp">cls</span><span class="p">,</span> <span class="n">s</span><span class="p">):</span>
        <span class="sa">r</span><span class="sd">&quot;&quot;&quot;</span>
<span class="sd">        Parse a Miller index string</span>

<span class="sd">        Negative indices can be denoted by:</span>
<span class="sd">         1. leading minus sign, e.g. ``[11-2]``</span>
<span class="sd">         2. trailing ``b`` (for &#39;bar&#39;), e.g. ``112b``</span>
<span class="sd">         3. LaTeX ``\bar{}``, e.g. ``[11\bar{2}]`` (which renders as :math:`[11\bar{2}]` in LaTeX)</span>

<span class="sd">        Leading or trailing brackets of various kinds are ignored.</span>
<span class="sd">        i.e. ``[001]``, ``{001}``, ``(001)``, ``[001]``, ``&lt;001&gt;``, ``001`` are all equivalent.</span>

<span class="sd">        Returns an array of components (i,j,k) or (h,k,i,l)</span>
<span class="sd">        &quot;&quot;&quot;</span>

        <span class="k">if</span> <span class="ow">not</span> <span class="nb">isinstance</span><span class="p">(</span><span class="n">s</span><span class="p">,</span> <span class="n">basestring</span><span class="p">):</span>
            <span class="k">raise</span> <span class="ne">TypeError</span><span class="p">(</span><span class="s2">&quot;Can&#39;t parse from </span><span class="si">%r</span><span class="s2"> of type </span><span class="si">%r</span><span class="s2">&quot;</span> <span class="o">%</span> <span class="p">(</span><span class="n">s</span><span class="p">,</span> <span class="nb">type</span><span class="p">(</span><span class="n">s</span><span class="p">)))</span>

        <span class="n">orig_s</span> <span class="o">=</span> <span class="n">s</span>
        <span class="k">for</span> <span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">b</span><span class="p">)</span> <span class="ow">in</span> <span class="p">[(</span><span class="sa">r</span><span class="s1">&#39;\bar{&#39;</span><span class="p">,</span><span class="s1">&#39;-&#39;</span><span class="p">)]</span> <span class="o">+</span> <span class="p">[(</span><span class="n">b</span><span class="p">,</span><span class="s1">&#39;&#39;</span><span class="p">)</span> <span class="k">for</span> <span class="n">b</span> <span class="ow">in</span> <span class="n">MillerIndex</span><span class="o">.</span><span class="n">all_brackets</span><span class="p">]:</span>
            <span class="n">s</span> <span class="o">=</span> <span class="n">s</span><span class="o">.</span><span class="n">replace</span><span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">b</span><span class="p">)</span>

        <span class="n">L</span> <span class="o">=</span> <span class="nb">list</span><span class="p">(</span><span class="n">s</span><span class="p">)</span>
        <span class="n">components</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">array</span><span class="p">([</span><span class="mi">1</span><span class="p">,</span><span class="mi">1</span><span class="p">,</span><span class="mi">1</span><span class="p">,</span><span class="mi">1</span><span class="p">])</span> <span class="c1"># space for up to 4 elements</span>
        <span class="n">i</span> <span class="o">=</span> <span class="mi">3</span> <span class="c1"># parse backwards from end of string</span>
        <span class="k">while</span> <span class="n">L</span><span class="p">:</span>
            <span class="k">if</span> <span class="n">i</span> <span class="o">&lt;</span> <span class="mi">0</span><span class="p">:</span>
                <span class="k">raise</span> <span class="ne">ValueError</span><span class="p">(</span><span class="s1">&#39;Cannot parse Miller index from string &quot;</span><span class="si">%s</span><span class="s1">&quot;, too many components found&#39;</span> <span class="o">%</span> <span class="n">orig_s</span><span class="p">)</span>
            <span class="n">c</span> <span class="o">=</span> <span class="n">L</span><span class="o">.</span><span class="n">pop</span><span class="p">()</span>
            <span class="k">if</span> <span class="n">c</span> <span class="o">==</span> <span class="s1">&#39;-&#39;</span><span class="p">:</span>
                <span class="k">if</span> <span class="n">i</span> <span class="o">==</span> <span class="mi">3</span><span class="p">:</span>
                    <span class="k">raise</span> <span class="ne">ValueError</span><span class="p">(</span><span class="s1">&#39;Miller index string &quot;</span><span class="si">%s</span><span class="s1">&quot; cannot end with a minus sign&#39;</span> <span class="o">%</span> <span class="n">orig_s</span><span class="p">)</span>
                <span class="n">components</span><span class="p">[</span><span class="n">i</span><span class="o">+</span><span class="mi">1</span><span class="p">]</span> <span class="o">*=</span> <span class="o">-</span><span class="mi">1</span>
            <span class="k">elif</span> <span class="n">c</span> <span class="o">==</span> <span class="s1">&#39;b&#39;</span><span class="p">:</span>
                <span class="n">components</span><span class="p">[</span><span class="n">i</span><span class="p">]</span> <span class="o">*=</span> <span class="o">-</span><span class="mi">1</span>
            <span class="k">elif</span> <span class="n">c</span> <span class="ow">in</span> <span class="p">[</span><span class="s1">&#39;0&#39;</span><span class="p">,</span> <span class="s1">&#39;1&#39;</span><span class="p">,</span> <span class="s1">&#39;2&#39;</span><span class="p">,</span> <span class="s1">&#39;3&#39;</span><span class="p">,</span> <span class="s1">&#39;4&#39;</span><span class="p">,</span> <span class="s1">&#39;5&#39;</span><span class="p">,</span> <span class="s1">&#39;6&#39;</span><span class="p">,</span> <span class="s1">&#39;7&#39;</span><span class="p">,</span> <span class="s1">&#39;8&#39;</span><span class="p">,</span> <span class="s1">&#39;9&#39;</span><span class="p">]:</span>
                <span class="n">components</span><span class="p">[</span><span class="n">i</span><span class="p">]</span> <span class="o">*=</span> <span class="nb">int</span><span class="p">(</span><span class="n">c</span><span class="p">)</span>
                <span class="n">i</span> <span class="o">-=</span> <span class="mi">1</span>
            <span class="k">else</span><span class="p">:</span>
                <span class="k">raise</span> <span class="ne">ValueError</span><span class="p">(</span><span class="s1">&#39;Unexpected character &quot;</span><span class="si">%s</span><span class="s1">&quot; in miller index string &quot;</span><span class="si">%s</span><span class="s1">&quot;&#39;</span> <span class="o">%</span> <span class="p">(</span><span class="n">c</span><span class="p">,</span> <span class="n">orig_s</span><span class="p">))</span>
                
        <span class="k">if</span> <span class="n">i</span> <span class="o">==</span> <span class="mi">0</span><span class="p">:</span>
            <span class="k">return</span> <span class="n">components</span><span class="p">[</span><span class="mi">1</span><span class="p">:]</span>
        <span class="k">elif</span> <span class="n">i</span> <span class="o">==</span> <span class="o">-</span><span class="mi">1</span><span class="p">:</span>
            <span class="k">return</span> <span class="n">components</span>
        <span class="k">else</span><span class="p">:</span>
            <span class="k">raise</span> <span class="ne">ValueError</span><span class="p">(</span><span class="s1">&#39;Cannot parse Miller index from string </span><span class="si">%s</span><span class="s1">, too few components found&#39;</span> <span class="o">%</span> <span class="n">orig_s</span><span class="p">)</span>
        
        <span class="bp">self</span><span class="o">.</span><span class="n">simplify</span><span class="p">()</span></div>

<div class="viewcode-block" id="MillerIndex.simplify"><a class="viewcode-back" href="../../structures.html#quippy.structures.MillerIndex.simplify">[docs]</a>    <span class="k">def</span> <span class="nf">simplify</span><span class="p">(</span><span class="bp">self</span><span class="p">):</span>
        <span class="sd">&quot;&quot;&quot;</span>
<span class="sd">        Simplify by dividing through by greatest common denominator</span>
<span class="sd">        &quot;&quot;&quot;</span>
        <span class="n">d</span> <span class="o">=</span> <span class="nb">abs</span><span class="p">(</span><span class="n">reduce</span><span class="p">(</span><span class="n">gcd</span><span class="p">,</span> <span class="bp">self</span><span class="p">))</span>
        <span class="bp">self</span><span class="p">[:]</span> <span class="o">/=</span> <span class="n">d</span></div>

    <span class="k">def</span> <span class="nf">simplified</span><span class="p">(</span><span class="bp">self</span><span class="p">):</span>
        <span class="n">copy</span> <span class="o">=</span> <span class="bp">self</span><span class="o">.</span><span class="n">copy</span><span class="p">()</span>
        <span class="n">copy</span><span class="o">.</span><span class="n">simplify</span><span class="p">()</span>
        <span class="k">return</span> <span class="n">copy</span>

<div class="viewcode-block" id="MillerIndex.normalised"><a class="viewcode-back" href="../../structures.html#quippy.structures.MillerIndex.normalised">[docs]</a>    <span class="k">def</span> <span class="nf">normalised</span><span class="p">(</span><span class="bp">self</span><span class="p">):</span>
        <span class="n">a</span> <span class="o">=</span> <span class="bp">self</span><span class="o">.</span><span class="n">as3</span><span class="p">()</span>
        <span class="k">return</span> <span class="n">a</span><span class="o">.</span><span class="n">copy</span><span class="p">()</span><span class="o">.</span><span class="n">view</span><span class="p">(</span><span class="n">FortranArray</span><span class="p">)</span><span class="o">/</span><span class="n">a</span><span class="o">.</span><span class="n">norm</span><span class="p">()</span></div>

    <span class="n">hat</span> <span class="o">=</span> <span class="n">normalised</span>

    <span class="k">def</span> <span class="nf">cross</span><span class="p">(</span><span class="bp">self</span><span class="p">,</span> <span class="n">other</span><span class="p">):</span>
        <span class="n">a</span> <span class="o">=</span> <span class="bp">self</span><span class="o">.</span><span class="n">as3</span><span class="p">()</span>
        <span class="n">b</span> <span class="o">=</span> <span class="n">MillerIndex</span><span class="p">(</span><span class="n">other</span><span class="p">)</span><span class="o">.</span><span class="n">as3</span><span class="p">()</span>
        <span class="k">return</span> <span class="n">np</span><span class="o">.</span><span class="n">cross</span><span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">b</span><span class="p">)</span><span class="o">.</span><span class="n">view</span><span class="p">(</span><span class="n">MillerIndex</span><span class="p">)</span><span class="o">.</span><span class="n">simplified</span><span class="p">()</span>

    <span class="k">def</span> <span class="nf">cosine</span><span class="p">(</span><span class="bp">self</span><span class="p">,</span> <span class="n">other</span><span class="p">):</span>
        <span class="n">other</span> <span class="o">=</span> <span class="n">MillerIndex</span><span class="p">(</span><span class="n">other</span><span class="p">)</span>
        <span class="k">return</span> <span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="bp">self</span><span class="o">.</span><span class="n">normalised</span><span class="p">(),</span> <span class="n">other</span><span class="o">.</span><span class="n">normalised</span><span class="p">())</span>

    <span class="k">def</span> <span class="nf">angle</span><span class="p">(</span><span class="bp">self</span><span class="p">,</span> <span class="n">other</span><span class="p">):</span>
        <span class="k">return</span> <span class="n">np</span><span class="o">.</span><span class="n">arccos</span><span class="p">(</span><span class="bp">self</span><span class="o">.</span><span class="n">cosine</span><span class="p">(</span><span class="n">other</span><span class="p">))</span>

    <span class="k">def</span> <span class="nf">as4</span><span class="p">(</span><span class="bp">self</span><span class="p">):</span>
        <span class="k">if</span> <span class="nb">len</span><span class="p">(</span><span class="bp">self</span><span class="p">)</span> <span class="o">==</span> <span class="mi">4</span><span class="p">:</span>
            <span class="k">return</span> <span class="bp">self</span>
        <span class="k">else</span><span class="p">:</span>
            <span class="n">h</span><span class="p">,</span> <span class="n">k</span><span class="p">,</span> <span class="n">l</span> <span class="o">=</span> <span class="bp">self</span>
            <span class="n">i</span> <span class="o">=</span> <span class="o">-</span><span class="p">(</span><span class="n">h</span><span class="o">+</span><span class="n">l</span><span class="p">)</span>
            <span class="k">return</span> <span class="n">MillerIndex</span><span class="p">((</span><span class="n">h</span><span class="p">,</span><span class="n">k</span><span class="p">,</span><span class="n">i</span><span class="p">,</span><span class="n">l</span><span class="p">))</span>

    <span class="k">def</span> <span class="nf">as3</span><span class="p">(</span><span class="bp">self</span><span class="p">):</span>
        <span class="k">if</span> <span class="nb">len</span><span class="p">(</span><span class="bp">self</span><span class="p">)</span> <span class="o">==</span> <span class="mi">3</span><span class="p">:</span>
            <span class="k">return</span> <span class="bp">self</span>
        <span class="k">else</span><span class="p">:</span>
            <span class="n">h</span><span class="p">,</span> <span class="n">k</span><span class="p">,</span> <span class="n">i</span><span class="p">,</span> <span class="n">l</span> <span class="o">=</span> <span class="bp">self</span>
            <span class="k">return</span> <span class="n">MillerIndex</span><span class="p">((</span><span class="n">h</span><span class="p">,</span> <span class="n">k</span><span class="p">,</span> <span class="n">l</span><span class="p">))</span>

    <span class="k">def</span> <span class="nf">plane_spacing</span><span class="p">(</span><span class="bp">self</span><span class="p">,</span> <span class="n">a</span><span class="p">):</span>
        <span class="k">return</span> <span class="n">a</span><span class="o">/</span><span class="bp">self</span><span class="o">.</span><span class="n">as3</span><span class="p">()</span><span class="o">.</span><span class="n">norm</span><span class="p">()</span></div>

<div class="viewcode-block" id="MillerPlane"><a class="viewcode-back" href="../../structures.html#quippy.structures.MillerPlane">[docs]</a><span class="k">def</span> <span class="nf">MillerPlane</span><span class="p">(</span><span class="n">v</span><span class="p">):</span>
   <span class="sd">&quot;&quot;&quot;Special case of :class:`MillerIndex` with ``type=&quot;plane&quot;``&quot;&quot;&quot;</span>
   <span class="k">return</span> <span class="n">MillerIndex</span><span class="p">(</span><span class="n">v</span><span class="p">,</span> <span class="s1">&#39;plane&#39;</span><span class="p">)</span></div>

<div class="viewcode-block" id="MillerDirection"><a class="viewcode-back" href="../../structures.html#quippy.structures.MillerDirection">[docs]</a><span class="k">def</span> <span class="nf">MillerDirection</span><span class="p">(</span><span class="n">v</span><span class="p">):</span>
   <span class="sd">&quot;&quot;&quot;Special case of :class:`MillerIndex` with ``type=&quot;direction&quot;`` (the default)&quot;&quot;&quot;</span>
   <span class="k">return</span> <span class="n">MillerIndex</span><span class="p">(</span><span class="n">v</span><span class="p">,</span> <span class="s1">&#39;direction&#39;</span><span class="p">)</span></div>


<div class="viewcode-block" id="angle_between"><a class="viewcode-back" href="../../structures.html#quippy.structures.angle_between">[docs]</a><span class="k">def</span> <span class="nf">angle_between</span><span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">b</span><span class="p">):</span>
    <span class="sd">&quot;&quot;&quot;Angle between crystallographic directions between a=[ijk] and b=[lmn], in radians.&quot;&quot;&quot;</span>
    <span class="k">return</span> <span class="n">MillerIndex</span><span class="p">(</span><span class="n">a</span><span class="p">)</span><span class="o">.</span><span class="n">angle</span><span class="p">(</span><span class="n">b</span><span class="p">)</span></div>


<div class="viewcode-block" id="rotation_matrix"><a class="viewcode-back" href="../../structures.html#quippy.structures.rotation_matrix">[docs]</a><span class="k">def</span> <span class="nf">rotation_matrix</span><span class="p">(</span><span class="n">unit</span><span class="p">,</span> <span class="n">y</span><span class="p">,</span> <span class="n">z</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span> <span class="n">x</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span> <span class="n">tol</span><span class="o">=</span><span class="mf">1e-5</span><span class="p">):</span>
    <span class="sd">&quot;&quot;&quot;Return 3x3 matrix rotation matrix defining a crack with open</span>
<span class="sd">    surface defined by the plane `y`=(l,m.n) or (h,k,i,l), and either</span>
<span class="sd">    crack tip line `z` or crack propagation direction `x`.&quot;&quot;&quot;</span>

    <span class="n">axes</span> <span class="o">=</span> <span class="n">fzeros</span><span class="p">((</span><span class="mi">3</span><span class="p">,</span><span class="mi">3</span><span class="p">))</span>
    <span class="n">y</span> <span class="o">=</span> <span class="n">MillerIndex</span><span class="p">(</span><span class="n">y</span><span class="p">)</span><span class="o">.</span><span class="n">as3</span><span class="p">()</span>

    <span class="k">if</span> <span class="p">(</span><span class="n">x</span> <span class="ow">is</span> <span class="kc">None</span> <span class="ow">and</span> <span class="n">z</span> <span class="ow">is</span> <span class="kc">None</span><span class="p">)</span> <span class="ow">or</span> <span class="p">(</span><span class="n">x</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span> <span class="ow">and</span> <span class="n">z</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span><span class="p">):</span>
        <span class="k">raise</span> <span class="ne">ValueError</span><span class="p">(</span><span class="s1">&#39;exactly one of x and z must be non-null&#39;</span><span class="p">)</span>

    <span class="n">axes</span><span class="p">[:,</span><span class="mi">2</span><span class="p">]</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="n">unit</span><span class="o">.</span><span class="n">g</span><span class="o">.</span><span class="n">T</span><span class="p">,</span> <span class="n">y</span><span class="p">)</span>     <span class="c1"># plane defined by y=(lmn)</span>

    <span class="k">if</span> <span class="n">z</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span><span class="p">:</span>
        <span class="n">z</span> <span class="o">=</span> <span class="n">MillerIndex</span><span class="p">(</span><span class="n">z</span><span class="p">)</span><span class="o">.</span><span class="n">as3</span><span class="p">()</span>
        <span class="n">axes</span><span class="p">[:,</span><span class="mi">3</span><span class="p">]</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="n">unit</span><span class="o">.</span><span class="n">lattice</span><span class="p">,</span> <span class="n">z</span><span class="p">)</span> <span class="c1"># line defined by z=[pqr]</span>

        <span class="n">axes</span><span class="p">[:,</span><span class="mi">2</span><span class="p">]</span> <span class="o">=</span> <span class="n">axes</span><span class="p">[:,</span><span class="mi">2</span><span class="p">]</span><span class="o">/</span><span class="n">axes</span><span class="p">[:,</span><span class="mi">2</span><span class="p">]</span><span class="o">.</span><span class="n">norm</span><span class="p">()</span>
        <span class="n">axes</span><span class="p">[:,</span><span class="mi">3</span><span class="p">]</span> <span class="o">=</span> <span class="n">axes</span><span class="p">[:,</span><span class="mi">3</span><span class="p">]</span><span class="o">/</span><span class="n">axes</span><span class="p">[:,</span><span class="mi">3</span><span class="p">]</span><span class="o">.</span><span class="n">norm</span><span class="p">()</span>

        <span class="k">if</span> <span class="nb">abs</span><span class="p">(</span><span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="n">axes</span><span class="p">[:,</span><span class="mi">2</span><span class="p">],</span> <span class="n">axes</span><span class="p">[:,</span><span class="mi">3</span><span class="p">]))</span> <span class="o">&gt;</span> <span class="n">tol</span><span class="p">:</span>
            <span class="k">raise</span> <span class="ne">ValueError</span><span class="p">(</span><span class="s1">&#39;y (</span><span class="si">%s</span><span class="s1">) and z (</span><span class="si">%s</span><span class="s1">) directions are not perpendicular&#39;</span> <span class="o">%</span> <span class="p">(</span><span class="n">y</span><span class="p">,</span><span class="n">z</span><span class="p">))</span>

        <span class="n">axes</span><span class="p">[:,</span><span class="mi">1</span><span class="p">]</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">cross</span><span class="p">(</span><span class="n">axes</span><span class="p">[:,</span><span class="mi">2</span><span class="p">],</span> <span class="n">axes</span><span class="p">[:,</span><span class="mi">3</span><span class="p">])</span>
    <span class="k">else</span><span class="p">:</span>
        <span class="n">x</span> <span class="o">=</span> <span class="n">MillerIndex</span><span class="p">(</span><span class="n">x</span><span class="p">)</span><span class="o">.</span><span class="n">as3</span><span class="p">()</span>
        <span class="n">axes</span><span class="p">[:,</span><span class="mi">1</span><span class="p">]</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="n">unit</span><span class="o">.</span><span class="n">lattice</span><span class="p">,</span> <span class="n">x</span><span class="p">)</span>

        <span class="n">axes</span><span class="p">[:,</span><span class="mi">2</span><span class="p">]</span> <span class="o">=</span> <span class="n">axes</span><span class="p">[:,</span><span class="mi">2</span><span class="p">]</span><span class="o">/</span><span class="n">axes</span><span class="p">[:,</span><span class="mi">2</span><span class="p">]</span><span class="o">.</span><span class="n">norm</span><span class="p">()</span>
        <span class="n">axes</span><span class="p">[:,</span><span class="mi">1</span><span class="p">]</span> <span class="o">=</span> <span class="n">axes</span><span class="p">[:,</span><span class="mi">1</span><span class="p">]</span><span class="o">/</span><span class="n">axes</span><span class="p">[:,</span><span class="mi">1</span><span class="p">]</span><span class="o">.</span><span class="n">norm</span><span class="p">()</span>

        <span class="k">if</span> <span class="nb">abs</span><span class="p">(</span><span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="n">axes</span><span class="p">[:,</span><span class="mi">2</span><span class="p">],</span> <span class="n">axes</span><span class="p">[:,</span><span class="mi">3</span><span class="p">]))</span> <span class="o">&gt;</span> <span class="n">tol</span><span class="p">:</span>
            <span class="k">raise</span> <span class="ne">ValueError</span><span class="p">(</span><span class="s1">&#39;y (</span><span class="si">%s</span><span class="s1">) and x (</span><span class="si">%s</span><span class="s1">) directions are not perpendicular&#39;</span> <span class="o">%</span> <span class="p">(</span><span class="n">y</span><span class="p">,</span><span class="n">x</span><span class="p">))</span>

        <span class="n">axes</span><span class="p">[:,</span><span class="mi">3</span><span class="p">]</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">cross</span><span class="p">(</span><span class="n">axes</span><span class="p">[:,</span><span class="mi">1</span><span class="p">],</span> <span class="n">axes</span><span class="p">[:,</span><span class="mi">2</span><span class="p">])</span>

    <span class="c1"># Rotation matrix is transpose of axes matrix</span>
    <span class="k">return</span> <span class="n">axes</span><span class="o">.</span><span class="n">T</span></div>



<div class="viewcode-block" id="orthorhombic_slab"><a class="viewcode-back" href="../../structures.html#quippy.structures.orthorhombic_slab">[docs]</a><span class="k">def</span> <span class="nf">orthorhombic_slab</span><span class="p">(</span><span class="n">at</span><span class="p">,</span> <span class="n">tol</span><span class="o">=</span><span class="mf">1e-5</span><span class="p">,</span> <span class="n">min_nrep</span><span class="o">=</span><span class="mi">1</span><span class="p">,</span> <span class="n">max_nrep</span><span class="o">=</span><span class="mi">5</span><span class="p">,</span> <span class="n">graphics</span><span class="o">=</span><span class="kc">False</span><span class="p">,</span> <span class="n">rot</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span> <span class="n">periodicity</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span> <span class="n">vacuum</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span> <span class="n">shift</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span> <span class="n">verbose</span><span class="o">=</span><span class="kc">True</span><span class="p">):</span>
    <span class="sd">&quot;&quot;&quot;Try to construct an orthorhombic cell equivalent to the</span>
<span class="sd">       primitive cell `at`, using supercells up to at most `max_nrep`</span>
<span class="sd">       repeats. Symmetry must be exact within a tolerance of `tol`. If</span>
<span class="sd">       `rot` is not None, we first transform `at` by the rotation</span>
<span class="sd">       matrix `rot`. The optional argument `periodicity` can be used to</span>
<span class="sd">       fix the periodicity one or more directions. It should be a three</span>
<span class="sd">       component vector with value zero in the unconstrained</span>
<span class="sd">       directions. The vector `vacuum` can be used to add vacuum in one</span>
<span class="sd">       or more directions. `shift` is a three component vector which</span>
<span class="sd">       can be used to shift the positions in the final cell. &quot;&quot;&quot;</span>

    <span class="k">def</span> <span class="nf">atoms_near_plane</span><span class="p">(</span><span class="n">at</span><span class="p">,</span> <span class="n">n</span><span class="p">,</span> <span class="n">d</span><span class="p">,</span> <span class="n">tol</span><span class="o">=</span><span class="mf">1e-5</span><span class="p">):</span>
        <span class="sd">&quot;&quot;&quot;Return a list of atoms within a distance `tol` of the plane defined by np.dot(n, at.pos) == d&quot;&quot;&quot;</span>
        <span class="n">pd</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="n">n</span><span class="p">,</span> <span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">)</span> <span class="o">-</span> <span class="n">d</span>
        <span class="k">return</span> <span class="p">(</span><span class="nb">abs</span><span class="p">(</span><span class="n">pd</span><span class="p">)</span> <span class="o">&lt;</span> <span class="n">tol</span><span class="p">)</span><span class="o">.</span><span class="n">nonzero</span><span class="p">()[</span><span class="mi">0</span><span class="p">]</span>

    <span class="k">def</span> <span class="nf">sort_by_distance</span><span class="p">(</span><span class="n">at</span><span class="p">,</span> <span class="n">ref_atom</span><span class="p">,</span> <span class="nb">dir</span><span class="p">,</span> <span class="n">candidates</span><span class="p">):</span>
        <span class="sd">&quot;&quot;&quot;Return a copy of `candidates` sorted by perpendicular distance from `ref_atom` in direction `dir`&quot;&quot;&quot;</span>
        <span class="n">distances_candidates</span> <span class="o">=</span>  <span class="nb">zip</span><span class="p">([</span><span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">i</span><span class="p">]</span><span class="o">-</span><span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">ref_atom</span><span class="p">]</span> <span class="k">for</span> <span class="n">i</span> <span class="ow">in</span> <span class="n">candidates</span><span class="p">],</span> <span class="n">candidates</span><span class="p">)</span>
        <span class="n">distances_candidates</span><span class="o">.</span><span class="n">sort</span><span class="p">()</span>
        <span class="k">return</span> <span class="p">[</span><span class="n">p</span> <span class="k">for</span> <span class="p">(</span><span class="n">d</span><span class="p">,</span> <span class="n">p</span><span class="p">)</span> <span class="ow">in</span> <span class="n">distances_candidates</span><span class="p">]</span>

    <span class="k">def</span> <span class="nf">orthorhombic_box</span><span class="p">(</span><span class="n">at</span><span class="p">):</span>
        <span class="sd">&quot;&quot;&quot;Return a copy of `at` in an orthorhombic box surrounded by vacuum&quot;&quot;&quot;</span>
        <span class="n">at</span> <span class="o">=</span> <span class="n">at</span><span class="o">.</span><span class="n">copy</span><span class="p">()</span>
        <span class="n">at</span><span class="o">.</span><span class="n">map_into_cell</span><span class="p">()</span>
        <span class="n">at</span><span class="o">.</span><span class="n">set_lattice</span><span class="p">([[</span><span class="mf">2.0</span><span class="o">*</span><span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="mi">1</span><span class="p">,:]</span><span class="o">.</span><span class="n">max</span><span class="p">()</span> <span class="o">-</span> <span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="mi">1</span><span class="p">,:]</span><span class="o">.</span><span class="n">min</span><span class="p">()),</span> <span class="mf">0.0</span><span class="p">,</span> <span class="mf">0.0</span><span class="p">],</span>
                        <span class="p">[</span><span class="mf">0.0</span><span class="p">,</span> <span class="mf">2.0</span><span class="o">*</span><span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="mi">2</span><span class="p">,:]</span><span class="o">.</span><span class="n">max</span><span class="p">()</span> <span class="o">-</span> <span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="mi">2</span><span class="p">,:]</span><span class="o">.</span><span class="n">min</span><span class="p">()),</span> <span class="mf">0.0</span><span class="p">],</span>
                        <span class="p">[</span><span class="mf">0.0</span><span class="p">,</span> <span class="mf">0.0</span><span class="p">,</span> <span class="mf">2.0</span><span class="o">*</span><span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="mi">3</span><span class="p">,:]</span><span class="o">.</span><span class="n">max</span><span class="p">()</span> <span class="o">-</span> <span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="mi">3</span><span class="p">,:]</span><span class="o">.</span><span class="n">min</span><span class="p">())]],</span>
                       <span class="n">scale_positions</span><span class="o">=</span><span class="kc">False</span><span class="p">)</span>
        <span class="n">at</span><span class="o">.</span><span class="n">map_into_cell</span><span class="p">()</span>
        <span class="k">return</span> <span class="n">at</span>

    <span class="k">def</span> <span class="nf">discard_outliers</span><span class="p">(</span><span class="n">at</span><span class="p">,</span> <span class="n">indices</span><span class="p">,</span> <span class="n">dirs</span><span class="p">,</span> <span class="n">keep_fraction</span><span class="o">=</span><span class="mf">0.5</span><span class="p">):</span>
        <span class="sd">&quot;&quot;&quot;Return a copy of `indices` with the atoms with fractional coordinates along directions in `dirs`</span>
<span class="sd">           outside +/-`keep_fraction`/2 excluded. Lattice used is close fitting, `at.lattice`/2.&quot;&quot;&quot;</span>
        <span class="n">g</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">linalg</span><span class="o">.</span><span class="n">inv</span><span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">lattice</span><span class="o">/</span><span class="mi">2</span><span class="p">)</span>
        <span class="n">t</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="n">g</span><span class="p">,</span> <span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">[:,</span><span class="n">indices</span><span class="p">])</span>
        <span class="k">return</span> <span class="n">indices</span><span class="p">[</span> <span class="n">np</span><span class="o">.</span><span class="n">logical_and</span><span class="p">(</span><span class="n">t</span><span class="p">[</span><span class="n">dirs</span><span class="p">,:]</span> <span class="o">&gt;=</span> <span class="o">-</span><span class="n">keep_fraction</span><span class="o">/</span><span class="mf">2.0</span><span class="p">,</span> <span class="n">t</span><span class="p">[</span><span class="n">dirs</span><span class="p">,:]</span> <span class="o">&lt;</span> <span class="n">keep_fraction</span><span class="o">/</span><span class="mf">2.0</span><span class="p">)</span><span class="o">.</span><span class="n">all</span><span class="p">(</span><span class="n">axis</span><span class="o">=</span><span class="mi">1</span><span class="p">)</span> <span class="p">]</span>

    <span class="k">def</span> <span class="nf">check_candidate_plane</span><span class="p">(</span><span class="n">at</span><span class="p">,</span> <span class="n">ref_plane</span><span class="p">,</span> <span class="n">cand_plane</span><span class="p">,</span> <span class="n">dirs</span><span class="p">,</span> <span class="n">verbose</span><span class="o">=</span><span class="kc">False</span><span class="p">,</span> <span class="n">label</span><span class="o">=</span><span class="s1">&#39;&#39;</span><span class="p">):</span>
        <span class="sd">&quot;&quot;&quot;Check whether in-plane displacements of atoms listed in `ref_plane` match those of `cand_plane` in directions given by `dirs`&quot;&quot;&quot;</span>

        <span class="c1"># Which pair of planes has more atoms, reference or candidate?</span>
        <span class="k">if</span> <span class="nb">len</span><span class="p">(</span><span class="n">ref_plane</span><span class="p">)</span> <span class="o">&lt;</span> <span class="nb">len</span><span class="p">(</span><span class="n">cand_plane</span><span class="p">):</span>
            <span class="n">smaller</span> <span class="o">=</span> <span class="n">ref_plane</span>
            <span class="n">larger</span>  <span class="o">=</span> <span class="n">cand_plane</span>
        <span class="k">else</span><span class="p">:</span>
            <span class="n">smaller</span> <span class="o">=</span> <span class="n">cand_plane</span>
            <span class="n">larger</span>  <span class="o">=</span> <span class="n">ref_plane</span>

        <span class="n">matches</span> <span class="o">=</span> <span class="p">{}</span>
        <span class="k">for</span> <span class="n">j</span> <span class="ow">in</span> <span class="n">smaller</span><span class="p">:</span>
            <span class="k">for</span> <span class="n">k</span> <span class="ow">in</span> <span class="n">larger</span><span class="p">:</span>
                <span class="k">if</span> <span class="n">at</span><span class="o">.</span><span class="n">z</span><span class="p">[</span><span class="n">k</span><span class="p">]</span> <span class="o">==</span> <span class="n">at</span><span class="o">.</span><span class="n">z</span><span class="p">[</span><span class="n">j</span><span class="p">]</span> <span class="ow">and</span> <span class="nb">abs</span><span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="n">dirs</span><span class="p">,</span><span class="n">k</span><span class="p">]</span> <span class="o">-</span> <span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="n">dirs</span><span class="p">,</span><span class="n">j</span><span class="p">])</span><span class="o">.</span><span class="n">max</span><span class="p">()</span> <span class="o">&lt;</span> <span class="n">tol</span><span class="p">:</span>
                    <span class="n">matches</span><span class="p">[</span><span class="n">j</span><span class="p">]</span> <span class="o">=</span> <span class="n">k</span>
                    <span class="k">break</span>

        <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
            <span class="nb">print</span> <span class="s1">&#39;   &#39;</span><span class="p">,</span> <span class="n">label</span><span class="p">,</span> <span class="nb">len</span><span class="p">(</span><span class="n">matches</span><span class="p">),</span> <span class="s1">&#39;/&#39;</span><span class="p">,</span> <span class="nb">len</span><span class="p">(</span><span class="n">smaller</span><span class="p">),</span> <span class="s1">&#39;matches&#39;</span>

        <span class="k">return</span> <span class="nb">len</span><span class="p">(</span><span class="n">matches</span><span class="p">)</span> <span class="o">==</span> <span class="nb">len</span><span class="p">(</span><span class="n">smaller</span><span class="p">)</span>


    <span class="k">if</span> <span class="n">rot</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span><span class="p">:</span>
        <span class="n">at</span> <span class="o">=</span> <span class="n">transform</span><span class="p">(</span><span class="n">at</span><span class="p">,</span> <span class="n">rot</span><span class="p">)</span>

    <span class="n">xyz</span> <span class="o">=</span> <span class="n">fidentity</span><span class="p">(</span><span class="mi">3</span><span class="p">)</span>
    <span class="n">nrep</span> <span class="o">=</span> <span class="n">min_nrep</span><span class="o">-</span><span class="mi">1</span>
    <span class="n">max_dist</span> <span class="o">=</span> <span class="n">fzeros</span><span class="p">(</span><span class="mi">3</span><span class="p">)</span>

    <span class="k">if</span> <span class="n">periodicity</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span><span class="p">:</span>
        <span class="n">periodicity</span> <span class="o">=</span> <span class="n">farray</span><span class="p">(</span><span class="n">periodicity</span><span class="p">)</span>
        <span class="n">periodicity</span> <span class="o">=</span> <span class="nb">dict</span><span class="p">(</span><span class="nb">zip</span><span class="p">((</span><span class="n">periodicity</span> <span class="o">!=</span> <span class="mi">0</span><span class="p">)</span><span class="o">.</span><span class="n">nonzero</span><span class="p">()[</span><span class="mi">0</span><span class="p">],</span> <span class="n">periodicity</span><span class="p">[</span><span class="n">periodicity</span> <span class="o">!=</span> <span class="mi">0</span><span class="p">]))</span>
    <span class="k">else</span><span class="p">:</span>
        <span class="n">periodicity</span> <span class="o">=</span> <span class="p">{}</span>

    <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
        <span class="k">for</span> <span class="p">(</span><span class="nb">dir</span><span class="p">,</span> <span class="n">p</span><span class="p">)</span> <span class="ow">in</span> <span class="n">periodicity</span><span class="o">.</span><span class="n">iteritems</span><span class="p">():</span>
            <span class="nb">print</span> <span class="s1">&#39;Periodicity in direction </span><span class="si">%d</span><span class="s1"> fixed at </span><span class="si">%f</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="p">(</span><span class="nb">dir</span><span class="p">,</span> <span class="n">p</span><span class="p">)</span>

    <span class="k">if</span> <span class="n">graphics</span><span class="p">:</span>
        <span class="kn">import</span> <span class="nn">atomeye</span>
        <span class="n">viewer</span> <span class="o">=</span> <span class="n">atomeye</span><span class="o">.</span><span class="n">AtomEyeViewer</span><span class="p">()</span>

    <span class="k">while</span> <span class="nb">sorted</span><span class="p">(</span><span class="n">periodicity</span><span class="o">.</span><span class="n">keys</span><span class="p">())</span> <span class="o">!=</span> <span class="p">[</span><span class="mi">1</span><span class="p">,</span><span class="mi">2</span><span class="p">,</span><span class="mi">3</span><span class="p">]:</span>
        <span class="n">nrep</span> <span class="o">+=</span> <span class="mi">1</span>
        <span class="k">if</span> <span class="n">nrep</span> <span class="o">&gt;</span> <span class="n">max_nrep</span><span class="p">:</span>
            <span class="k">raise</span> <span class="ne">ValueError</span><span class="p">(</span><span class="s1">&#39;Maximum size of supercell (</span><span class="si">%d</span><span class="s1">) exceeded&#39;</span> <span class="o">%</span> <span class="n">max_nrep</span><span class="p">)</span>
        <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
            <span class="nb">print</span> <span class="s1">&#39;</span><span class="se">\n\n</span><span class="s1">Supercell </span><span class="si">%d</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="n">nrep</span>
        <span class="n">sup</span> <span class="o">=</span> <span class="n">supercell</span><span class="p">(</span><span class="n">at</span><span class="p">,</span> <span class="n">nrep</span><span class="p">,</span> <span class="n">nrep</span><span class="p">,</span> <span class="n">nrep</span><span class="p">)</span>
        <span class="n">box</span> <span class="o">=</span> <span class="n">orthorhombic_box</span><span class="p">(</span><span class="n">sup</span><span class="p">)</span>
        <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[:]</span> <span class="o">=</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span> <span class="o">-</span> <span class="n">np</span><span class="o">.</span><span class="n">tile</span><span class="p">(</span><span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="o">.</span><span class="n">mean</span><span class="p">(</span><span class="n">axis</span><span class="o">=</span><span class="mi">2</span><span class="p">),</span> <span class="p">[</span><span class="n">box</span><span class="o">.</span><span class="n">n</span><span class="p">,</span> <span class="mi">1</span><span class="p">])</span><span class="o">.</span><span class="n">T</span>

        <span class="k">for</span> <span class="nb">dir</span> <span class="ow">in</span> <span class="nb">set</span><span class="p">([</span><span class="mi">1</span><span class="p">,</span><span class="mi">2</span><span class="p">,</span><span class="mi">3</span><span class="p">])</span> <span class="o">-</span> <span class="nb">set</span><span class="p">(</span><span class="n">periodicity</span><span class="o">.</span><span class="n">keys</span><span class="p">()):</span>

            <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
                <span class="nb">print</span> <span class="s1">&#39;  Direction </span><span class="si">%d</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="nb">dir</span>

            <span class="n">other_dirs</span> <span class="o">=</span> <span class="nb">list</span><span class="p">(</span><span class="nb">set</span><span class="p">([</span><span class="mi">1</span><span class="p">,</span><span class="mi">2</span><span class="p">,</span><span class="mi">3</span><span class="p">])</span> <span class="o">-</span> <span class="nb">set</span><span class="p">([</span><span class="nb">dir</span><span class="p">]))</span>

            <span class="n">pos_index</span> <span class="o">=</span> <span class="nb">zip</span><span class="p">(</span><span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,:],</span><span class="n">frange</span><span class="p">(</span><span class="n">box</span><span class="o">.</span><span class="n">n</span><span class="p">))</span>
            <span class="n">pos_index</span><span class="o">.</span><span class="n">sort</span><span class="p">()</span>

            <span class="c1"># Find a pair of planes</span>
            <span class="k">while</span> <span class="n">pos_index</span><span class="p">:</span>
                <span class="n">ref_pos1</span><span class="p">,</span> <span class="n">ref_atom1</span> <span class="o">=</span> <span class="n">pos_index</span><span class="o">.</span><span class="n">pop</span><span class="p">(</span><span class="mi">0</span><span class="p">)</span>

                <span class="c1"># Find atom to define second plane</span>
                <span class="k">while</span> <span class="n">pos_index</span><span class="p">:</span>
                    <span class="n">ref_pos2</span><span class="p">,</span> <span class="n">ref_atom2</span> <span class="o">=</span> <span class="n">pos_index</span><span class="o">.</span><span class="n">pop</span><span class="p">(</span><span class="mi">0</span><span class="p">)</span>
                    <span class="k">if</span> <span class="nb">abs</span><span class="p">(</span><span class="n">ref_pos2</span> <span class="o">-</span> <span class="n">ref_pos1</span><span class="p">)</span> <span class="o">&gt;</span> <span class="n">tol</span><span class="p">:</span> <span class="k">break</span>
                <span class="k">else</span><span class="p">:</span>
                    <span class="k">continue</span>

                <span class="n">ref_plane1</span> <span class="o">=</span> <span class="n">atoms_near_plane</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">xyz</span><span class="p">[:,</span><span class="nb">dir</span><span class="p">],</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">ref_atom1</span><span class="p">],</span> <span class="n">tol</span><span class="p">)</span>
                <span class="n">ref_plane2</span> <span class="o">=</span> <span class="n">atoms_near_plane</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">xyz</span><span class="p">[:,</span><span class="nb">dir</span><span class="p">],</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">ref_atom2</span><span class="p">],</span> <span class="n">tol</span><span class="p">)</span>

                <span class="c1"># Only keep reference atoms in the centre of the cell</span>
                <span class="n">ref_plane1</span> <span class="o">=</span> <span class="n">discard_outliers</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">ref_plane1</span><span class="p">,</span> <span class="n">other_dirs</span><span class="p">)</span>
                <span class="n">ref_plane2</span> <span class="o">=</span> <span class="n">discard_outliers</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">ref_plane2</span><span class="p">,</span> <span class="n">other_dirs</span><span class="p">)</span>

                <span class="k">if</span> <span class="nb">len</span><span class="p">(</span><span class="n">ref_plane1</span><span class="p">)</span> <span class="o">&gt;</span> <span class="mi">2</span> <span class="ow">and</span> <span class="nb">len</span><span class="p">(</span><span class="n">ref_plane2</span><span class="p">)</span> <span class="o">&gt;</span> <span class="mi">2</span><span class="p">:</span>
                    <span class="c1"># Now we&#39;ve got two planes, both with more than two atoms in them</span>
                    <span class="k">break</span>
            <span class="k">else</span><span class="p">:</span>
                <span class="c1"># Used up all atoms without finding two planes</span>
                <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
                    <span class="nb">print</span> <span class="s1">&#39;    No valid reference planes found.</span><span class="se">\n</span><span class="s1">&#39;</span>
                <span class="k">continue</span>

            <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
                <span class="nb">print</span> <span class="s1">&#39;    Reference plane #1 through atom </span><span class="si">%d</span><span class="s1"> &#39;</span> <span class="o">%</span> <span class="n">ref_atom1</span>
                <span class="nb">print</span> <span class="s1">&#39;    Reference plane #2 through atom </span><span class="si">%d</span><span class="s1"> at distance </span><span class="si">%r</span><span class="se">\n</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="p">(</span><span class="n">ref_atom2</span><span class="p">,</span> <span class="n">ref_pos2</span> <span class="o">-</span> <span class="n">ref_pos1</span><span class="p">)</span>

            <span class="k">if</span> <span class="n">graphics</span><span class="p">:</span>
                <span class="n">highlight</span> <span class="o">=</span> <span class="n">fzeros</span><span class="p">(</span><span class="n">box</span><span class="o">.</span><span class="n">n</span><span class="p">)</span>
                <span class="n">highlight</span><span class="p">[</span><span class="n">ref_plane1</span><span class="p">]</span> <span class="o">=</span> <span class="mi">1</span>
                <span class="n">highlight</span><span class="p">[</span><span class="n">ref_plane2</span><span class="p">]</span> <span class="o">=</span> <span class="mi">2</span>
                <span class="n">box</span><span class="o">.</span><span class="n">add_property</span><span class="p">(</span><span class="s1">&#39;highlight&#39;</span><span class="p">,</span> <span class="n">highlight</span><span class="p">,</span> <span class="n">overwrite</span><span class="o">=</span><span class="kc">True</span><span class="p">)</span>
                <span class="n">viewer</span><span class="o">.</span><span class="n">show</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="s1">&#39;highlight&#39;</span><span class="p">)</span>
                <span class="n">viewer</span><span class="o">.</span><span class="n">wait</span><span class="p">()</span>
                <span class="n">raw_input</span><span class="p">(</span><span class="s1">&#39;continue&#39;</span><span class="p">)</span>

            <span class="n">candidates</span> <span class="o">=</span> <span class="p">[</span><span class="n">i</span> <span class="k">for</span> <span class="n">i</span> <span class="ow">in</span> <span class="n">frange</span><span class="p">(</span><span class="n">box</span><span class="o">.</span><span class="n">n</span><span class="p">)</span> <span class="k">if</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">i</span><span class="p">]</span> <span class="o">&gt;</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">ref_atom2</span><span class="p">]</span> <span class="o">+</span> <span class="n">max_dist</span><span class="p">[</span><span class="nb">dir</span><span class="p">]</span> <span class="o">+</span> <span class="n">tol</span><span class="p">]</span>
            <span class="n">candidates</span> <span class="o">=</span> <span class="n">sort_by_distance</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">ref_atom1</span><span class="p">,</span> <span class="nb">dir</span><span class="p">,</span> <span class="n">candidates</span><span class="p">)</span>

            <span class="k">while</span> <span class="n">candidates</span><span class="p">:</span>
                <span class="n">cand1</span> <span class="o">=</span> <span class="n">candidates</span><span class="o">.</span><span class="n">pop</span><span class="p">(</span><span class="mi">0</span><span class="p">)</span>

                <span class="n">max_dist</span><span class="p">[</span><span class="nb">dir</span><span class="p">]</span> <span class="o">=</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">cand1</span><span class="p">]</span> <span class="o">-</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">ref_atom1</span><span class="p">]</span>

                <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
                    <span class="nb">print</span> <span class="s1">&#39;    Trying plane through atom </span><span class="si">%d</span><span class="s1"> distance </span><span class="si">%r</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="p">(</span><span class="n">cand1</span><span class="p">,</span> <span class="n">max_dist</span><span class="p">[</span><span class="nb">dir</span><span class="p">])</span>

                <span class="n">cand_plane1</span> <span class="o">=</span> <span class="n">atoms_near_plane</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">xyz</span><span class="p">[:,</span><span class="nb">dir</span><span class="p">],</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">cand1</span><span class="p">],</span> <span class="n">tol</span><span class="p">)</span>

                <span class="k">for</span> <span class="n">cand2</span> <span class="ow">in</span> <span class="n">sort_by_distance</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">ref_atom1</span><span class="p">,</span> <span class="nb">dir</span><span class="p">,</span> <span class="nb">set</span><span class="p">(</span><span class="n">candidates</span><span class="p">)</span> <span class="o">-</span> <span class="nb">set</span><span class="p">(</span><span class="n">cand_plane1</span><span class="p">)):</span>
                    <span class="k">if</span> <span class="nb">abs</span><span class="p">((</span><span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">cand2</span><span class="p">]</span> <span class="o">-</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">cand1</span><span class="p">])</span> <span class="o">-</span> <span class="p">(</span><span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">ref_atom2</span><span class="p">]</span> <span class="o">-</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">ref_atom1</span><span class="p">]))</span> <span class="o">&lt;</span> <span class="n">tol</span><span class="p">:</span>
                        <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
                            <span class="nb">print</span> <span class="s1">&#39;    Found pair to plane, passing through atom </span><span class="si">%d</span><span class="s1"> distance </span><span class="si">%r</span><span class="s1"> &#39;</span> <span class="o">%</span> <span class="p">(</span><span class="n">cand2</span><span class="p">,</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">cand2</span><span class="p">]</span> <span class="o">-</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">ref_atom1</span><span class="p">])</span>
                        <span class="k">break</span>
                <span class="k">else</span><span class="p">:</span>
                    <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
                        <span class="nb">print</span> <span class="s1">&#39;    Cannot find second candidate plane.</span><span class="se">\n</span><span class="s1">&#39;</span>
                    <span class="n">candidates</span> <span class="o">=</span> <span class="n">sort_by_distance</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">ref_atom1</span><span class="p">,</span> <span class="nb">dir</span><span class="p">,</span> <span class="nb">set</span><span class="p">(</span><span class="n">candidates</span><span class="p">)</span> <span class="o">-</span> <span class="nb">set</span><span class="p">(</span><span class="n">cand_plane1</span><span class="p">))</span>
                    <span class="k">continue</span>

                <span class="k">if</span> <span class="n">graphics</span><span class="p">:</span>
                    <span class="n">highlight</span><span class="p">[</span><span class="n">cand_plane1</span><span class="p">]</span> <span class="o">=</span> <span class="mi">3</span>
                    <span class="n">box</span><span class="o">.</span><span class="n">highlight</span><span class="p">[:]</span> <span class="o">=</span> <span class="n">highlight</span>
                    <span class="n">viewer</span><span class="o">.</span><span class="n">show</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="s1">&#39;highlight&#39;</span><span class="p">)</span>
                    <span class="n">viewer</span><span class="o">.</span><span class="n">wait</span><span class="p">()</span>

                <span class="n">cand_plane2</span> <span class="o">=</span> <span class="n">atoms_near_plane</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">xyz</span><span class="p">[:,</span><span class="nb">dir</span><span class="p">],</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">cand2</span><span class="p">],</span> <span class="n">tol</span><span class="p">)</span>

                <span class="k">if</span> <span class="n">graphics</span><span class="p">:</span>
                    <span class="n">highlight</span><span class="p">[</span><span class="n">cand_plane2</span><span class="p">]</span> <span class="o">=</span> <span class="mi">4</span>
                    <span class="n">box</span><span class="o">.</span><span class="n">highlight</span><span class="p">[:]</span> <span class="o">=</span> <span class="n">highlight</span>                    
                    <span class="n">viewer</span><span class="o">.</span><span class="n">show</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="s1">&#39;highlight&#39;</span><span class="p">)</span>
                    <span class="n">viewer</span><span class="o">.</span><span class="n">wait</span><span class="p">()</span>

                    <span class="n">highlight</span><span class="p">[</span><span class="n">cand_plane1</span><span class="p">]</span> <span class="o">=</span> <span class="mi">0</span>
                    <span class="n">highlight</span><span class="p">[</span><span class="n">cand_plane2</span><span class="p">]</span> <span class="o">=</span> <span class="mi">0</span>

                <span class="c1"># Remove cand_plane1 from list of candidates</span>
                <span class="n">candidates</span> <span class="o">=</span> <span class="n">sort_by_distance</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">ref_atom1</span><span class="p">,</span> <span class="nb">dir</span><span class="p">,</span> <span class="nb">set</span><span class="p">(</span><span class="n">candidates</span><span class="p">)</span> <span class="o">-</span> <span class="nb">set</span><span class="p">(</span><span class="n">cand_plane1</span><span class="p">))</span>

                <span class="c1"># Check ref_plane1 against cand_plane1 and ref_plane2 against cand_plane2 in directions</span>
                <span class="c1"># listed in other_dirs</span>
                <span class="n">match1</span> <span class="o">=</span> <span class="n">check_candidate_plane</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">ref_plane1</span><span class="p">,</span> <span class="n">cand_plane1</span><span class="p">,</span> <span class="n">other_dirs</span><span class="p">,</span> <span class="n">verbose</span><span class="p">,</span> <span class="s1">&#39;Plane #1:&#39;</span><span class="p">)</span>
                <span class="n">match2</span> <span class="o">=</span> <span class="n">check_candidate_plane</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="n">ref_plane2</span><span class="p">,</span> <span class="n">cand_plane2</span><span class="p">,</span> <span class="n">other_dirs</span><span class="p">,</span> <span class="n">verbose</span><span class="p">,</span> <span class="s1">&#39;Plane #2:&#39;</span><span class="p">)</span>

                <span class="k">if</span> <span class="n">match1</span> <span class="ow">and</span> <span class="n">match2</span><span class="p">:</span>
                    <span class="n">periodicity</span><span class="p">[</span><span class="nb">dir</span><span class="p">]</span> <span class="o">=</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">cand1</span><span class="p">]</span> <span class="o">-</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">ref_atom1</span><span class="p">]</span>
                    <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
                        <span class="nb">print</span> <span class="s1">&#39;</span><span class="se">\n</span><span class="s1">  Periodicity in direction </span><span class="si">%d</span><span class="s1"> is </span><span class="si">%f</span><span class="se">\n</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="p">(</span><span class="nb">dir</span><span class="p">,</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">cand1</span><span class="p">]</span> <span class="o">-</span> <span class="n">box</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="nb">dir</span><span class="p">,</span><span class="n">ref_atom1</span><span class="p">])</span>

                    <span class="k">if</span> <span class="n">graphics</span><span class="p">:</span>
                        <span class="n">highlight</span><span class="p">[</span><span class="n">cand_plane1</span><span class="p">]</span> <span class="o">=</span> <span class="mi">3</span>
                        <span class="n">highlight</span><span class="p">[</span><span class="n">cand_plane2</span><span class="p">]</span> <span class="o">=</span> <span class="mi">3</span>
                        <span class="n">box</span><span class="o">.</span><span class="n">highlight</span><span class="p">[:]</span> <span class="o">=</span> <span class="n">highlight</span>
                        <span class="n">viewer</span><span class="o">.</span><span class="n">show</span><span class="p">(</span><span class="n">box</span><span class="p">,</span> <span class="s1">&#39;highlight&#39;</span><span class="p">)</span>
                        <span class="n">viewer</span><span class="o">.</span><span class="n">wait</span><span class="p">()</span>
                        <span class="n">raw_input</span><span class="p">(</span><span class="s1">&#39;continue...&#39;</span><span class="p">)</span>
                    <span class="k">break</span>

                <span class="k">if</span> <span class="n">graphics</span><span class="p">:</span>
                    <span class="n">raw_input</span><span class="p">(</span><span class="s1">&#39;continue...&#39;</span><span class="p">)</span>
            <span class="k">else</span><span class="p">:</span>
                <span class="c1"># Failed to find match for direction dir</span>
                <span class="k">continue</span>

    <span class="c1"># Finally, construct new cell by selecting atoms within first unit cell</span>
    <span class="n">lattice</span> <span class="o">=</span> <span class="n">farray</span><span class="p">(</span><span class="n">np</span><span class="o">.</span><span class="n">diag</span><span class="p">([</span><span class="n">periodicity</span><span class="p">[</span><span class="mi">1</span><span class="p">],</span> <span class="n">periodicity</span><span class="p">[</span><span class="mi">2</span><span class="p">],</span> <span class="n">periodicity</span><span class="p">[</span><span class="mi">3</span><span class="p">]]))</span>
    <span class="n">g</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">linalg</span><span class="o">.</span><span class="n">inv</span><span class="p">(</span><span class="n">lattice</span><span class="p">)</span>

    <span class="n">nrepx</span><span class="p">,</span> <span class="n">nrepy</span><span class="p">,</span> <span class="n">nrepz</span> <span class="o">=</span> <span class="n">fit_box_in_cell</span><span class="p">(</span><span class="n">periodicity</span><span class="p">[</span><span class="mi">1</span><span class="p">],</span> <span class="n">periodicity</span><span class="p">[</span><span class="mi">2</span><span class="p">],</span> <span class="n">periodicity</span><span class="p">[</span><span class="mi">3</span><span class="p">],</span> <span class="n">at</span><span class="o">.</span><span class="n">lattice</span><span class="p">)</span>

    <span class="n">sup</span> <span class="o">=</span> <span class="n">supercell</span><span class="p">(</span><span class="n">at</span><span class="p">,</span> <span class="n">nrepx</span><span class="p">,</span> <span class="n">nrepy</span><span class="p">,</span> <span class="n">nrepz</span><span class="p">)</span>
    <span class="n">sup</span><span class="o">.</span><span class="n">map_into_cell</span><span class="p">()</span>

    <span class="c1"># small shift to avoid coincidental cell alignments</span>
    <span class="n">delta</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">tile</span><span class="p">([</span><span class="mf">0.01</span><span class="p">,</span> <span class="mf">0.02</span><span class="p">,</span> <span class="mf">0.03</span><span class="p">],</span> <span class="p">[</span><span class="n">sup</span><span class="o">.</span><span class="n">n</span><span class="p">,</span> <span class="mi">1</span><span class="p">])</span><span class="o">.</span><span class="n">T</span>
    <span class="k">if</span> <span class="n">shift</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span> <span class="ow">and</span> <span class="n">vacuum</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span><span class="p">:</span>
        <span class="n">delta</span> <span class="o">=</span> <span class="n">delta</span> <span class="o">+</span> <span class="n">np</span><span class="o">.</span><span class="n">tile</span><span class="p">(</span><span class="n">shift</span><span class="p">,</span> <span class="p">[</span><span class="n">sup</span><span class="o">.</span><span class="n">n</span><span class="p">,</span> <span class="mi">1</span><span class="p">])</span><span class="o">.</span><span class="n">T</span>
    <span class="n">t</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="n">g</span><span class="p">,</span> <span class="n">sup</span><span class="o">.</span><span class="n">pos</span><span class="p">)</span> <span class="o">+</span> <span class="n">delta</span>

    <span class="n">orthorhombic</span> <span class="o">=</span> <span class="n">sup</span><span class="o">.</span><span class="n">select</span><span class="p">(</span><span class="n">np</span><span class="o">.</span><span class="n">logical_and</span><span class="p">(</span><span class="n">t</span> <span class="o">&gt;=</span> <span class="o">-</span><span class="mf">0.5</span><span class="p">,</span> <span class="n">t</span> <span class="o">&lt;</span> <span class="mf">0.5</span><span class="p">)</span><span class="o">.</span><span class="n">all</span><span class="p">(</span><span class="n">axis</span><span class="o">=</span><span class="mi">1</span><span class="p">))</span>

    <span class="k">if</span> <span class="n">vacuum</span><span class="p">:</span>
        <span class="n">lattice</span> <span class="o">=</span> <span class="n">farray</span><span class="p">(</span><span class="n">np</span><span class="o">.</span><span class="n">diag</span><span class="p">(</span><span class="n">lattice</span><span class="o">.</span><span class="n">diagonal</span><span class="p">()</span> <span class="o">+</span> <span class="n">vacuum</span><span class="p">))</span>

    <span class="k">if</span> <span class="n">shift</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span> <span class="ow">and</span> <span class="n">vacuum</span> <span class="ow">is</span> <span class="kc">None</span><span class="p">:</span>
        <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
            <span class="nb">print</span> <span class="s1">&#39;Shifting positions by </span><span class="si">%s</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="n">lattice</span><span class="p">,</span> <span class="n">shift</span><span class="p">)</span>
        <span class="n">orthorhombic</span><span class="o">.</span><span class="n">pos</span> <span class="o">+=</span> <span class="n">np</span><span class="o">.</span><span class="n">tile</span><span class="p">(</span><span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="n">lattice</span><span class="p">,</span> <span class="n">shift</span><span class="p">),</span> <span class="p">[</span><span class="n">orthorhombic</span><span class="o">.</span><span class="n">n</span><span class="p">,</span> <span class="mi">1</span><span class="p">])</span><span class="o">.</span><span class="n">T</span>

    <span class="n">orthorhombic</span><span class="o">.</span><span class="n">set_lattice</span><span class="p">(</span><span class="n">lattice</span><span class="p">,</span> <span class="n">scale_positions</span><span class="o">=</span><span class="kc">False</span><span class="p">)</span>
    <span class="n">orthorhombic</span><span class="o">.</span><span class="n">map_into_cell</span><span class="p">()</span>
    <span class="k">return</span> <span class="n">orthorhombic</span></div>



<span class="n">quartz_params</span> <span class="o">=</span> <span class="p">{</span><span class="s1">&#39;experiment&#39;</span><span class="p">:</span> <span class="p">{</span><span class="s1">&#39;a&#39;</span><span class="p">:</span> <span class="mf">4.9160</span><span class="p">,</span>
                                <span class="s1">&#39;c&#39;</span><span class="p">:</span> <span class="mf">5.4054</span><span class="p">,</span>
                                <span class="s1">&#39;u&#39;</span><span class="p">:</span> <span class="mf">0.4697</span><span class="p">,</span>
                                <span class="s1">&#39;x&#39;</span><span class="p">:</span> <span class="mf">0.4135</span><span class="p">,</span>
                                <span class="s1">&#39;y&#39;</span><span class="p">:</span> <span class="mf">0.2669</span><span class="p">,</span>
                                <span class="s1">&#39;z&#39;</span><span class="p">:</span> <span class="mf">0.1191</span><span class="p">},</span>
                 <span class="s1">&#39;CASTEP_LDA&#39;</span><span class="p">:</span> <span class="p">{</span><span class="s1">&#39;a&#39;</span><span class="p">:</span> <span class="mf">4.87009</span><span class="p">,</span>
                                <span class="s1">&#39;c&#39;</span><span class="p">:</span> <span class="mf">5.36255</span><span class="p">,</span>
                                <span class="s1">&#39;u&#39;</span><span class="p">:</span> <span class="mf">0.46699</span><span class="p">,</span>
                                <span class="s1">&#39;x&#39;</span><span class="p">:</span> <span class="mf">0.41289</span><span class="p">,</span>
                                <span class="s1">&#39;y&#39;</span><span class="p">:</span> <span class="mf">0.27198</span><span class="p">,</span>
                                <span class="s1">&#39;z&#39;</span><span class="p">:</span> <span class="mf">0.11588</span><span class="p">},</span>
                 <span class="s1">&#39;CASTEP_GGA&#39;</span><span class="p">:</span> <span class="p">{</span><span class="s1">&#39;a&#39;</span><span class="p">:</span> <span class="mf">5.02836</span><span class="p">,</span>
                                <span class="s1">&#39;c&#39;</span><span class="p">:</span> <span class="mf">5.51193</span><span class="p">,</span>
                                <span class="s1">&#39;u&#39;</span><span class="p">:</span> <span class="mf">0.48128</span><span class="p">,</span>
                                <span class="s1">&#39;x&#39;</span><span class="p">:</span> <span class="mf">0.41649</span><span class="p">,</span>
                                <span class="s1">&#39;y&#39;</span><span class="p">:</span> <span class="mf">0.24661</span><span class="p">,</span>
                                <span class="s1">&#39;z&#39;</span><span class="p">:</span> <span class="mf">0.13594</span><span class="p">},</span>
                 <span class="s1">&#39;ASAP_JRK&#39;</span><span class="p">:</span> <span class="p">{</span><span class="s1">&#39;a&#39;</span><span class="p">:</span> <span class="mf">4.8403809707320216</span><span class="p">,</span>
                              <span class="s1">&#39;c&#39;</span><span class="p">:</span> <span class="mf">5.3285240037002248</span><span class="p">,</span>
                              <span class="s1">&#39;u&#39;</span><span class="p">:</span> <span class="mf">0.46417561617105912</span><span class="p">,</span>
                              <span class="s1">&#39;x&#39;</span><span class="p">:</span> <span class="mf">0.41174271054205958</span><span class="p">,</span>
                              <span class="s1">&#39;y&#39;</span><span class="p">:</span> <span class="mf">0.27872745399831672</span><span class="p">,</span>
                              <span class="s1">&#39;z&#39;</span><span class="p">:</span> <span class="mf">0.10973603276909905</span><span class="p">}</span>
                 <span class="p">}</span>

<span class="k">def</span> <span class="nf">get_quartz_params</span><span class="p">(</span><span class="n">at</span><span class="p">):</span>

    <span class="k">assert</span> <span class="n">at</span><span class="o">.</span><span class="n">n</span> <span class="o">==</span> <span class="mi">9</span>
    <span class="k">assert</span> <span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">z</span> <span class="o">==</span> <span class="mi">14</span><span class="p">)</span><span class="o">.</span><span class="n">sum</span><span class="p">()</span> <span class="o">==</span> <span class="mi">3</span>
    <span class="k">assert</span> <span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">z</span> <span class="o">==</span> <span class="mi">8</span><span class="p">)</span><span class="o">.</span><span class="n">sum</span><span class="p">()</span> <span class="o">==</span> <span class="mi">6</span>

    <span class="kn">from</span> <span class="nn">quippy</span> <span class="k">import</span> <span class="n">get_lattice_params</span>

    <span class="n">lat_params</span> <span class="o">=</span> <span class="n">get_lattice_params</span><span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">lattice</span><span class="p">)</span>
    <span class="n">a</span><span class="p">,</span> <span class="n">c</span> <span class="o">=</span> <span class="n">lat_params</span><span class="p">[</span><span class="mi">0</span><span class="p">],</span> <span class="n">lat_params</span><span class="p">[</span><span class="mi">2</span><span class="p">]</span>
    <span class="nb">print</span> <span class="s1">&#39;a      = &#39;</span><span class="p">,</span> <span class="n">a</span>
    <span class="nb">print</span> <span class="s1">&#39;c      = &#39;</span><span class="p">,</span> <span class="n">c</span>
    <span class="nb">print</span> <span class="s1">&#39;c/a    = &#39;</span><span class="p">,</span> <span class="n">c</span><span class="o">/</span><span class="n">a</span>
    <span class="nb">print</span> <span class="s1">&#39;V      = &#39;</span><span class="p">,</span> <span class="n">at</span><span class="o">.</span><span class="n">cell_volume</span><span class="p">()</span>
    <span class="nb">print</span> <span class="s1">&#39;V/SiO2 = &#39;</span><span class="p">,</span> <span class="n">at</span><span class="o">.</span><span class="n">cell_volume</span><span class="p">()</span><span class="o">/</span><span class="mf">3.0</span>

    <span class="n">frac_pos</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">dot</span><span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">g</span><span class="p">,</span> <span class="n">at</span><span class="o">.</span><span class="n">pos</span><span class="p">)</span>
    <span class="n">u</span> <span class="o">=</span> <span class="n">frac_pos</span><span class="p">[</span><span class="mi">1</span><span class="p">,</span><span class="mi">1</span><span class="p">]</span>
    <span class="n">x</span><span class="p">,</span><span class="n">y</span><span class="p">,</span><span class="n">z</span> <span class="o">=</span> <span class="n">frac_pos</span><span class="p">[:,</span><span class="mi">4</span><span class="p">]</span>
    <span class="n">z</span> <span class="o">-=</span> <span class="mf">2.0</span><span class="o">/</span><span class="mf">3.0</span>
    <span class="k">if</span> <span class="n">z</span> <span class="o">&lt;</span> <span class="mf">0.0</span><span class="p">:</span> <span class="n">z</span> <span class="o">+=</span> <span class="mf">1.0</span>
    <span class="k">if</span> <span class="n">z</span> <span class="o">&gt;</span> <span class="mf">1.0</span><span class="p">:</span> <span class="n">z</span> <span class="o">--</span> <span class="mf">1.0</span>

    <span class="nb">print</span> <span class="s1">&#39;u      = &#39;</span><span class="p">,</span> <span class="n">u</span>
    <span class="nb">print</span> <span class="s1">&#39;x      = &#39;</span><span class="p">,</span> <span class="n">x</span>
    <span class="nb">print</span> <span class="s1">&#39;y      = &#39;</span><span class="p">,</span> <span class="n">y</span>
    <span class="nb">print</span> <span class="s1">&#39;z      = &#39;</span><span class="p">,</span> <span class="n">z</span>

    <span class="k">return</span> <span class="p">{</span><span class="s1">&#39;a&#39;</span><span class="p">:</span><span class="n">a</span><span class="p">,</span> <span class="s1">&#39;c&#39;</span><span class="p">:</span><span class="n">c</span><span class="p">,</span> <span class="s1">&#39;u&#39;</span><span class="p">:</span><span class="n">u</span><span class="p">,</span> <span class="s1">&#39;x&#39;</span><span class="p">:</span><span class="n">x</span><span class="p">,</span> <span class="s1">&#39;y&#39;</span><span class="p">:</span><span class="n">y</span><span class="p">,</span> <span class="s1">&#39;z&#39;</span><span class="p">:</span><span class="n">z</span><span class="p">}</span>


<div class="viewcode-block" id="get_bond_lengths"><a class="viewcode-back" href="../../structures.html#quippy.structures.get_bond_lengths">[docs]</a><span class="k">def</span> <span class="nf">get_bond_lengths</span><span class="p">(</span><span class="n">at</span><span class="p">):</span>
    <span class="sd">&quot;&quot;&quot;Return a dictionary mapping tuples (Species1, Species2) to an farray of bond-lengths&quot;&quot;&quot;</span>
    <span class="n">at</span><span class="o">.</span><span class="n">calc_connect</span><span class="p">()</span>
    <span class="n">r_ij</span> <span class="o">=</span> <span class="n">farray</span><span class="p">(</span><span class="mf">0.0</span><span class="p">)</span>
    <span class="n">res</span> <span class="o">=</span> <span class="p">{}</span>
    <span class="k">for</span> <span class="n">i</span> <span class="ow">in</span> <span class="n">frange</span><span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">n</span><span class="p">):</span>
        <span class="k">for</span> <span class="n">n</span> <span class="ow">in</span> <span class="n">frange</span><span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">n_neighbours</span><span class="p">(</span><span class="n">i</span><span class="p">)):</span>
            <span class="n">j</span> <span class="o">=</span> <span class="n">at</span><span class="o">.</span><span class="n">neighbour</span><span class="p">(</span><span class="n">i</span><span class="p">,</span> <span class="n">n</span><span class="p">,</span> <span class="n">distance</span><span class="o">=</span><span class="n">r_ij</span><span class="p">)</span>
            <span class="nb">print</span> <span class="n">i</span><span class="p">,</span> <span class="n">j</span><span class="p">,</span> <span class="n">at</span><span class="o">.</span><span class="n">z</span><span class="p">[</span><span class="n">i</span><span class="p">],</span> <span class="n">at</span><span class="o">.</span><span class="n">z</span><span class="p">[</span><span class="n">j</span><span class="p">],</span> <span class="n">r_ij</span>
            <span class="n">minij</span><span class="p">,</span> <span class="n">maxij</span> <span class="o">=</span> <span class="nb">min</span><span class="p">((</span><span class="n">i</span><span class="p">,</span><span class="n">j</span><span class="p">)),</span> <span class="nb">max</span><span class="p">((</span><span class="n">i</span><span class="p">,</span><span class="n">j</span><span class="p">))</span>
            <span class="n">key</span> <span class="o">=</span> <span class="p">(</span><span class="n">at</span><span class="o">.</span><span class="n">species</span><span class="p">[</span><span class="n">minij</span><span class="p">]</span><span class="o">.</span><span class="n">stripstrings</span><span class="p">(),</span> <span class="n">at</span><span class="o">.</span><span class="n">species</span><span class="p">[</span><span class="n">maxij</span><span class="p">]</span><span class="o">.</span><span class="n">stripstrings</span><span class="p">())</span>
            <span class="k">if</span> <span class="ow">not</span> <span class="n">key</span> <span class="ow">in</span> <span class="n">res</span><span class="p">:</span> <span class="n">res</span><span class="p">[</span><span class="n">key</span><span class="p">]</span> <span class="o">=</span> <span class="p">[]</span>
            <span class="n">res</span><span class="p">[</span><span class="n">key</span><span class="p">]</span><span class="o">.</span><span class="n">append</span><span class="p">(</span><span class="n">r_ij</span><span class="o">.</span><span class="n">astype</span><span class="p">(</span><span class="nb">float</span><span class="p">))</span>
    <span class="nb">print</span> <span class="n">res</span>
    <span class="k">return</span> <span class="nb">dict</span><span class="p">((</span><span class="n">k</span><span class="p">,</span><span class="n">farray</span><span class="p">(</span><span class="n">v</span><span class="p">))</span> <span class="k">for</span> <span class="p">(</span><span class="n">k</span><span class="p">,</span><span class="n">v</span><span class="p">)</span> <span class="ow">in</span> <span class="n">res</span><span class="o">.</span><span class="n">iteritems</span><span class="p">())</span></div>

<div class="viewcode-block" id="get_bulk_params"><a class="viewcode-back" href="../../structures.html#quippy.structures.get_bulk_params">[docs]</a><span class="k">def</span> <span class="nf">get_bulk_params</span><span class="p">(</span><span class="n">bulk</span><span class="p">,</span> <span class="n">lat_type</span><span class="p">,</span> <span class="n">verbose</span><span class="o">=</span><span class="kc">True</span><span class="p">):</span>
    <span class="sd">&quot;&quot;&quot;Return 6-tuple of lattice parameters a, c, u, x, y, z for</span>
<span class="sd">       cell `bulk` of lattice type `lat_type`&quot;&quot;&quot;</span>
    <span class="n">a</span><span class="p">,</span> <span class="n">b</span><span class="p">,</span> <span class="n">c</span><span class="p">,</span> <span class="n">alpha</span><span class="p">,</span> <span class="n">beta</span><span class="p">,</span> <span class="n">gamma</span> <span class="o">=</span> <span class="n">get_lattice_params</span><span class="p">(</span><span class="n">bulk</span><span class="o">.</span><span class="n">lattice</span><span class="p">)</span>
    <span class="k">del</span> <span class="n">b</span><span class="p">,</span> <span class="n">alpha</span><span class="p">,</span> <span class="n">beta</span><span class="p">,</span> <span class="n">gamma</span>
    <span class="n">u</span><span class="p">,</span> <span class="n">x</span><span class="p">,</span> <span class="n">y</span><span class="p">,</span> <span class="n">z</span> <span class="o">=</span> <span class="p">(</span><span class="kc">None</span><span class="p">,</span> <span class="kc">None</span><span class="p">,</span> <span class="kc">None</span><span class="p">,</span> <span class="kc">None</span><span class="p">)</span>

    <span class="k">if</span> <span class="n">lat_type</span> <span class="ow">in</span> <span class="p">(</span><span class="s1">&#39;diamond&#39;</span><span class="p">,</span> <span class="s1">&#39;bcc&#39;</span><span class="p">,</span> <span class="s1">&#39;fcc&#39;</span><span class="p">):</span>
        <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
            <span class="nb">print</span> <span class="s1">&#39;</span><span class="si">%s</span><span class="s1"> lattice, a=</span><span class="si">%.3f</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="p">(</span><span class="n">lat_type</span><span class="p">,</span> <span class="n">a</span><span class="p">)</span>
    <span class="k">elif</span> <span class="n">lat_type</span> <span class="o">==</span> <span class="s1">&#39;anatase&#39;</span><span class="p">:</span>
        <span class="n">u</span> <span class="o">=</span> <span class="n">bulk</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="mi">3</span><span class="p">,</span> <span class="mi">5</span><span class="p">]</span><span class="o">/</span><span class="n">c</span>
        <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
            <span class="nb">print</span> <span class="s1">&#39;anatase lattice, a=</span><span class="si">%.3f</span><span class="s1"> c=</span><span class="si">%.3f</span><span class="s1"> u=</span><span class="si">%.3f</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">c</span><span class="p">,</span> <span class="n">u</span><span class="p">)</span>
    <span class="k">elif</span> <span class="n">lat_type</span> <span class="o">==</span> <span class="s1">&#39;rutile&#39;</span><span class="p">:</span>
        <span class="n">u</span> <span class="o">=</span> <span class="n">bulk</span><span class="o">.</span><span class="n">pos</span><span class="p">[</span><span class="mi">1</span><span class="p">,</span> <span class="mi">3</span><span class="p">]</span><span class="o">/</span><span class="n">a</span>
        <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
            <span class="nb">print</span> <span class="s1">&#39;rutile lattice, a=</span><span class="si">%.3f</span><span class="s1"> c=</span><span class="si">%.3f</span><span class="s1"> u=</span><span class="si">%.3f</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">c</span><span class="p">,</span> <span class="n">u</span><span class="p">)</span>
    <span class="k">elif</span> <span class="n">lat_type</span> <span class="o">==</span> <span class="s1">&#39;alpha_quartz&#39;</span><span class="p">:</span>
        <span class="n">qz</span> <span class="o">=</span> <span class="n">get_quartz_params</span><span class="p">(</span><span class="n">bulk</span><span class="p">)</span>
        <span class="n">a</span><span class="p">,</span> <span class="n">c</span><span class="p">,</span> <span class="n">u</span><span class="p">,</span> <span class="n">x</span><span class="p">,</span> <span class="n">y</span><span class="p">,</span> <span class="n">z</span> <span class="o">=</span> <span class="p">(</span><span class="n">qz</span><span class="p">[</span><span class="s1">&#39;a&#39;</span><span class="p">],</span> <span class="n">qz</span><span class="p">[</span><span class="s1">&#39;c&#39;</span><span class="p">],</span> <span class="n">qz</span><span class="p">[</span><span class="s1">&#39;u&#39;</span><span class="p">],</span>
                            <span class="n">qz</span><span class="p">[</span><span class="s1">&#39;x&#39;</span><span class="p">],</span> <span class="n">qz</span><span class="p">[</span><span class="s1">&#39;y&#39;</span><span class="p">],</span> <span class="n">qz</span><span class="p">[</span><span class="s1">&#39;z&#39;</span><span class="p">])</span>
        <span class="k">if</span> <span class="n">verbose</span><span class="p">:</span>
            <span class="nb">print</span> <span class="s1">&#39;alpha_quartz lattice, &#39;</span><span class="p">,</span>
            <span class="nb">print</span> <span class="p">(</span><span class="s1">&#39;a=</span><span class="si">%.3f</span><span class="s1"> c=</span><span class="si">%.3f</span><span class="s1"> u=</span><span class="si">%.3f</span><span class="s1"> x=</span><span class="si">%.3f</span><span class="s1"> y=</span><span class="si">%.3f</span><span class="s1"> z=</span><span class="si">%.3f</span><span class="s1">&#39;</span>
                   <span class="o">%</span> <span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">c</span><span class="p">,</span> <span class="n">u</span><span class="p">,</span> <span class="n">x</span><span class="p">,</span> <span class="n">y</span> <span class="p">,</span><span class="n">z</span><span class="p">))</span>
    <span class="k">else</span><span class="p">:</span>
        <span class="k">raise</span> <span class="ne">ValueError</span><span class="p">(</span><span class="s1">&#39;Unknown latttice type </span><span class="si">%s</span><span class="s1">&#39;</span> <span class="o">%</span> <span class="n">lat_type</span><span class="p">)</span>

    <span class="k">return</span> <span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">c</span><span class="p">,</span> <span class="n">u</span><span class="p">,</span> <span class="n">x</span><span class="p">,</span> <span class="n">y</span><span class="p">,</span> <span class="n">z</span><span class="p">)</span></div>
</pre></div>

           </div>
           
          </div>
          <footer>
  

  <hr/>

  <div role="contentinfo">
    <p>
        &copy; Copyright 2008-2016, James Kermode

    </p>
  </div>
  Built with <a href="http://sphinx-doc.org/">Sphinx</a> using a <a href="https://github.com/rtfd/sphinx_rtd_theme">theme</a> provided by <a href="https://readthedocs.org">Read the Docs</a>. 

</footer>

        </div>
      </div>

    </section>

  </div>
  


  <script type="text/javascript">
      jQuery(function () {
          SphinxRtdTheme.Navigation.enable(true);
      });
  </script>

  
  
    
   

</body>
</html>